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We investigate how the dynamics of a single chain influences the kinetics of early stage phase 
separation in a symmetric binary polymer mixture. We consider quenches from the disordered phase 
into the region of spinodal instability. On a mean field level we approach this problem with two 
j^jfj' methods: a dynamical extension of the self consistent field theory for Gaussian chains, with the 

r-{ i density variables evolving in time, and the method of the external potential dynamics where the 

■ effective external fields are propagated in time. Different wave vector dependencies of the kinetic 

coefficient are taken into account. These early stages of spinodal decomposition are also studied 
through Monte Carlo simulations employing the bond fluctuation model that maps the chains - 
in our case with 64 effective segments - on a coarse grained lattice. The results obtained through 
(— I | self consistent field calculations and Monte Carlo simulations can be compared because the time, 

. length, and temperature scales are mapped onto each other through the diffusion constant, the chain 

extension, and the energy of mixing. The quantitative comparison of the relaxation rate of the global 
structure factor shows that a kinetic coefficient according to the Rouse model gives a much better 
agreement than a local, i.e. wave vector independent, kinetic factor. Including fluctuations in the 
self consistent field calculations leads to a shorter time span of spinodal behaviour and a reduction of 
* , the relaxation rate for smaller wave vectors and prevents the relaxation rate from becoming negative 

for larger values of the wave vector. This is also in agreement with the simulation results. 
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I. INTRODUCTION 



Over many years the study of phase transitions in fluid mixtures has become an important field because of their 
omnipresent occurrence in nature and more importantly because of their technological meaning in creating new mate- 
rials [l]J^]: blending of different species can reduce cost, improve processibility, provide synergy between components 
OO I and allow for recycling. 

Sophisticated analytical techniques - e.g., self consistent field theory ^-^] or P-RISM theory Q - exist for calculating 
the phase behaviour and detailed interfacial properties of polymer blends. Because of the large extension of polymer 



chains equilibrium properties are well described; however, the analytical description is much less satisfactory for the 
dynamics of multi-component systems. 

The transition from a completely homogeneous mixture to an equilibrated two phase system is a very long process 
consisting of a sequence of highly inhomogeneous states. Many different methods to analyse this process have been 
applied including experiments, theory, and computational simulations. After quenching a polymer blend from the 
one phase region deep into the two phase region spinodal decomposition takes place. Different time regimes are 
recognised during spinodal decomposition. During early stages the amplitudes of the concentration fluctuations, that 
are amplified, are still so small that they do not interact with each other. For later times the local composition reaches 
the equilibrium value of the coexisting phases, and nonlinear interactions between the fluctuation modes become more 
£j , and more important. At even later stages hydrodynamics dominates the coarsening [fjj. 

Spinodal decomposition can be approximately described through time dependent Ginzburg-Landau theory, also 
known as Cahn-Hilliard-Cook theory ]9|-pd]| . Many efforts have been made to numerically calculate the time evolution 
of a mixture in the framework of this theory [p^|-pd[. These calculations are appropriate to give a first insight into 
how phase separation takes place, but are far from describing the system quantitatively well. Even though binary 
polymer blends are an ideal testing bed for these approaches, earlier Monte Carlo simulations jlB] [l7| found rather 
pronounced deviations from their predictions. To find a quantitatively better description, it is necessary to go beyond 
Ginzburg-Landau theory. Self consistent field theory (SCFT) has proven to be one of the most successful 

descriptions of equilibrium properties of polymer mixtures on a mean field level. The idea to use SCFT to develop a 
dynamical mean field theory is not new ]l8|-^l|], but, although SCFT describes equilibrium properties well, usually the 
influence of the single chain dynamics is neglected due to computational expenses, again leading to a more qualitative 
than quantitative description of the phase separation. In most calculations a simple constant Onsager coefficient 
is used which would be appropriate if the movement of the polymers were comparable to the movement of point 
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particles. Remembering the connectivity of polymers it is not likely that this local coupling is sufficient to lead to a 
quantitatively correct description of the dynamics and other nonlocal Onsager coefficients have been proposed [^2|~^4| . 

In this study we are interested in the early stages of spinodal decomposition in an incompressible symmetric binary 
polymer mixture after a quench from the one phase region to the two phase region paying special attention to the 
differences in the collective dynamics when the chains obey local and Rouse dynamics. We therefore employ two 
dynamical versions of self consistent field theory for polymer mixtures. The first method, referred to as dynamic self 
consistent field theory ]f8| , |l9t (DSCFT), uses the free energy functional of the SCFT for an incompressible mixture 
in terms of the local composition. This free energy functional leads to a kinetic equation for the time dependent 
local composition which is integrated numerically. This method has been applied for studying the ordering kinetics 
in block copolymer and surfactant systems JlJ,|25|]. We use both a local and a nonlocal Onsager transport coefficient, 
as is expected for polymers obeying Rouse dynamics. 

The underlying idea of the second method is to express the free energy of the investigated system only in terms of 
the effective external fields, which are thermodynamically conjugated to the composition. This description leads to 
an equation of motion for the effective external fields, which is integrated in time. Following Maurits et al. |2(| who 
derived a similar equation of motion neglecting random fluctuations we refer to this method as the external potential 
dynamics (EPD), which has been shown to automatically incorporate Rouse dynamics. 

The aim of our study is to explore both methods and to investigate the role of the Onsager coefficient. To 
decide whether our results resemble any 'real' dynamics a comparison with results obtained through other methods 
is preferable. Although many experiments have analysed spinodal decomposition in binary polymer mixtures only 
few [^7|-^9) have looked at early stages of demixing and the influence of the Onsager coefficient. Another difficulty in 
comparing experimental results with our calculations lies in the fact that actual demixing is influenced by many factors, 
for example, preparation of the probe, polydispersity, or strong dynamic asymmetries due to different glass transition 
temperatures of the two species ]3C|j3l| . Therefore a quantitative mapping between theory and experimental systems 
is difficult. To gain a quantitative insight we compare our mean field calculations with results obtained through Monte 
Carlo simulations. No parameters have to be adjusted for this comparison. In the simulations we employ the bond 
fluctuation model [j32] , |33| which is well established for studying properties of polymer melts. It has been shown that 
the equilibrium properties of a polymer mixture as given by the self consistent field theory are almost quantitatively 
reproduced through Monte Carlo simulations using the bond fluctuation model |34||35| . 

Our paper is organised in the following way. In section [ij we introduce the self consistent field theory for a binary 
polymer mixture, explaining how equilibrium properties can be derived with this method. Then we introduce general 
aspects of the dynamics in a polymer mixture and show how they can be incorporated in the SCFT leading to 
dynamic self consistent field theory. The following section presents another approach to including dynamics in SCFT 
regarding the effective external fields as the time dependent variable. Implications regarding density fluctuations 
in this description are discussed in detail. Section |l| serves as a brief introduction to the performed Monte Carlo 
simulations and explains how a direct comparison with SCFT calculations is possible. In section IV results of our 
calculations and simulations are presented first showing the difference between the use of local and nonlocal coupling, 
then comparing the mean field results with the simulations and finally analysing the influence of random fluctuations. 
The paper finishes with a summary. 



II. SELF CONSISTENT FIELD THEORY 



We consider an incompressible mixture of A and B polymers consisting of n polymers in a volume V = L x xL y xL z 
with periodic boundary conditions. There are ua polymers of kind A in the system, for B polymers respectively 
Hb = n — tia- In SCFT IHHm polymers are modelled as Gaussian chains with the end-to-end distance R e . In 
addition we choose the A and B polymers to both have the same number of monomers Na — Nb — N and to be of 
the same architecture. The overall particle density in the system is denoted as p = nN/V. The microscopic density 
cj>A of the A monomers can be defined through the polymer conformations {r Q (r)} : 

n nA r 1 

= - E / drS(r-r iA (r)) (I) 

< r < 1 parameterises the contour of a chain. For B monomers a similar equation holds. (All the following equations 
regarding only A monomers, are equivalently valid for B monomers without mentioning this explicitly.) Regarding a 
repulsion between the two kinds of polymers that is expressed through the Flory-Huggins parameter x, the canonical 
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partition function has the following form: 

Z rv. 



n A 
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The functional integral V sums over all possible conformations of the chains. Va[*] is the so called Wiener measure 



Va[y] ~ exp 



(3) 



which represents the statistical weight of a noninteracting Gaussian chain. Repulsive interactions at short interparticlc 
distances reduce fluctuations of the total density. This is incorporated effectively via the incompressibility constraint 
which is expressed through the (5-function in Eq.(|^) p7y . 

This defines a system with many interacting polymer chains; its partition function obviously cannot be analytically 
solved. Inserting new auxiliary field variables $>a, &b, Wa, Wb via a Hubbard-Stratonovich transformation, we can 
reformulate the many polymer problem in terms of a single polymer in external fields. 



Z ~ J V$ A VW A V$ B VW B 5($a + $b - 1) exp [-F [W A) W B ,<f> A , $ fl ] /k B T] 



(4) 



Thus we have found an expression for the canonical free energy depending on the new variables and the single chain 
partition function Q A : 



F[Wa,W b ,$a,®b] 
k B T 
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(5) 



<f> A = Nn A / pV denoting the average density of the A polymers in the system and 

fi 



V[r A ]V A [r] exp 



drW A {r{T)) 



(G) 



being the partition function of a single A chain in the external field Wa- 



A. Equilibrium 

The functional integral in Eq.(|j) also cannot be calculated explicitly. Therefore we employ a saddle point approx- 
imation. This means that only the largest contribution to the integrand is considered and the integration does not 
have to be carried out. The saddle point approximation of Eq.(^j) is equivalent to the minimisation of the free energy 
(|J) with respect to the auxiliary variables. The values of the fields and densities at the saddle point are denoted by 
lower case letters and are given by the set of equations that has to be solved self consistently: 
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(9) 
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= 6* A and 



Eqs.(^) and (|^) make clear that 
in equilibrium the average microscopic density ((f> A ) actually equals the thermal average monomer density of a single 



The averages of the microscopic densities are given by: 



vb) = <P B - 
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A chain in an external field wa- In other words after starting with a description where we would have to take all 
possible interactions of every polymer chain into account, see Eq.(|^), we have now found a mean field description of 
the system in which it is sufficient to look at a single chain in an effective external field. 

To calculate the monomer densities it is useful to define the end segment distribution (^(r, t), which gives the 
probability to find the end of a chain with length t at position r when exposed to a field wa'- 



- I dTw A (r(T)) 



ii 



(11) 



q A (r,t)= J V[r(t)]P A [r(t)]S(r(t)-T)exp 
For this end segment distribution the following diffusion equation holds []3j 

: R 2 e V 2 qA(r,t)-w A qA{r,t) (12 



dq A (r,t) _ 1 „ . 



dt 6 

with the boundary condition qA (r, 0) = 1. After finding a solution to this equation the monomer density is immediately 
given through the following expression : 

&to = ^ I dtq A (r,t)q A (T,l-t) (13) 
Qa Jo 

The single chain partition function can be calculated with: 



d 3 rq A {r,l) (14) 

v 



After having presented all necessary equations to calculate equilibrium properties of polymer mixtures, we are now 
interested in using this description for calculating the dynamics in a polymer mixture. Two ways to achieve this are 
introduced in the following sections. 



B. Dynamic Self Consistent Field Theory (DSCFT) 

Because of the fact that the concentrations of the polymers are conserved the continuity equation is valid: 

^> + VJ(r,t) = C> (15) 

J(r, t) denoting the current density of the monomers at position r at time t. One now assumes p2|-p4| a linear 
relation between the current density and the gradient of the exchange potential fi = — (Note there is only 
one independent chemical potential in the system because of the incompressibility constraint) 

J(r,t)=- [ dr' 3 A(r,r')V/i(r',i) (16) 
Jv 

with the kinetic coefficient A(r, r') describing the connection between the force acting on the monomers through the 
gradient of the chemical potential at position r' and the resulting current density at position r. This describes a 
purely relaxational dynamics, effects due to hydrodynamic flow are not captured. The Onsager coefficient A can be 
modelled in different ways: The simplest approach would be local coupling which results in the Onsager coefficient 
being proportional to the local density. Bearing in mind that polymers have a certain extension it is clear that 
nonlocal coupling should lead to a better description, although local coupling is often used in calculations of dynamic 
models based on Ginzburg-Landau type energy functionals for simplicity reasons p2|-pT| . In the Rouse model forces 



acting on a monomer caused by the other monomers are also taken into account |38| 
that is proportional to the pair-correlation function p2|- p4| , p6| . These two approac 
coefficients: 



39| ] . This leads to a kinetic factor 
les lead to the following Onsager 



Aiocai(r) = DNcj) A {r, t)<p B (r,t) local coupling (17) 

A Rousc (r,r') a DN$ a $ b P q (t,i') Rouse (18) 

D denoting the single chain diffusion constant and Po(r, r') the pair-correlation function. The Rouse Onsager coeffi- 
cient written in Eq. (|l8|) is only approximately valid, when the pair-correlation function is the same for both polymer 
species. A more general expression is found in reference Eq]. 
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Another model for nonlocal coupling is the reptation model |3^,Q which is appropriate for polymer melts with 
very long chains, i.e. very entangled chains. Hereby the idea is that a polymer chain is constrained by the other 
polymers and is forced to move along the polymer tube axis. This dynamics is expected p2|-p4[ to also lead to an 
Onsager coefficient that is proportional to the pair-correlation function. Therefore the influence of single chain Rouse 
and reptation dynamics on the collective dynamics of the system is qualitatively comparable. Reptation shall not be 
regarded in the further study. 

Eqs.([l5]) and ( |l6| ) lead to the following diffusion equation the last term representing noise that obeys the fluctuation- 
dissipation theorem. 



at 

After the Fourier transformation this diffusion equation has a simple form 

d(j} A (q,t) 



= V / dr ,3 A(r,r')V^(r',t) +r](r,t) (19) 
Jv 



at 



A(q)<rV(q,t)+7 7 (q,t) (20) 



These diffusion equations implicitly assume that the relaxation time of the chain conformations is smaller than the 
time scale on which composition fluctuations evolve. The chain conformations are expected to be 'in equilibrium' 
with respect to the instantaneous spatially varying composition. 

To use the introduced diffusion equation ([l9|) in the frame of SCFT, let us return to the general expression (||) that 
describes the canonical partition function depending upon the variables $a, Wa, Wb that are independent of each 
other. If we now employ a saddle point approximation in the variables Wa and Wb we find Eqs.(^) and (||), that de- 
scribe a unique relation between the fields wa, wb and the densities <pA, 4>b- This means we can calculate the densities 
explicitely with Eqs.(|l2]), (|l3|), and (|l4|), if the fields wa and Wb are known. Therefore with this approximation the free 
energy F becomes a function that depends on the densities (or fields) only: F [4>a,4>b,wa[4>a],wb[4>b\\ — F [4> Al 4>b\- 
With this free energy it is now possible to calculate the exchange potential /z(r): 

fi(r) 6F[<j> A {r),cl> B {r)] 6F[<p A {r),^ B (r)] 1 

= 72 i \ c ■ / x = T7 ix^Ww) ~ <Mr)) - [wa <M r ) ~ Wfl feW (21) 

k B T 54>A(r) 54> b {y) N 

The simple form of diffusion equation ( p0| ) and the diffusion equation for the end segment distribution of a polymer 
chain ( |i"2| ) suggests to use a Fourier expansion of all spatially dependent variables for actual calculations because the 
Fourier functions are eigenfunctions of the squared gradient V 2 . The following set of orthonormal functions is used 
in all our SCFT calculations: 



JlmnX?) = norm(t) norm(m) norm(n) cos x I cos y I cos z I I, m, n — 0, 1, 2, . 



L/<r: J \ Ly J \ L 7 



norm(«) = 



V2 i^O 
1 i = 



(22) 



We have now found all necessary equations to numerically calculate the time evolution of the densities in a binary 
polymer mixture, leading us the following procedure we refer to as the dynamic self consistent field theory (DSCFT) 
method: First we have a given density profile at time t — 0. As mentioned before Eqs.([l2]), ( |l3| ) and ( fl4| ) give us the 
possibility to explicitly calculate the single chain densities for known external fields. Unfortunately the inversion is 
needed because we have given densities and have to find the matching external fields. This leads to a set of nonlinear 
equations, that is numerically solved through the Newton-Broyden method [[0J. After the external fields, that 'create' 
the given density profile, have been found the exchange potential [i is calculated via Eq. ( pl| ) and then inserted into the 
diffusion equation (^0|). The diffusion equation is subsequently integrated using a simplified Runge-Kutta method. 
This leads to a new given density after a discrete time step and the whole procedure starts anew. 

Apart from regarding local dynamics in our DSCFT calculations it is favourable to consider nonlocal coupling 
because this leads to a better description of the dynamics in a polymer mixture. However in DSCFT the difficulty 
in using Rouse dynamics lies in the computational expense of calculating the pair-correlation function for each time 
step. As an approximation during early stages of dembdng the pair-correlation function of a homogeneous melt, as 
it is given through the random phase approximation (RPA) |^|, is used leading to the following Onsager coefficient: 

A(g) = DNM B 2{X + e ~ X ~ 1) (23) 
x z 

x is defined as x = i? 2 q 2 /6 with R e denoting the end-to-end distance of a polymer. 



5 



C. External Potential Dynamics 



In the previously introduced DSCFT method we needed a way to reduce the number of independent variables in the 
partition function. Through the saddle point approximation in the fields we obtained a free energy functional in terms 
of the densities. This, in turn, yields a Langcvin equation for the dynamics of the densities. For the external potential 
dynamics (EPD) method [ p6| we are looking for a way to express the dynamics of the binary polymer mixture through 
an equation of motion for the external fields W a and Wb ■ Our starting point is again the canonical partition function 
(^). Via a Hubbard-Stratonovich transformation we introduce the the field variables W = Wa~Wb and U = Wa+Wb 
and obtain: 



Z 



n A 
■ exp 



1 „ / n A n B \ 

y| HI! V[v ia \D[v 1b ]Va[v ia \Pb[v 1b ]\ 
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NJv 



d 6 v 



U 



This defines a free energy function G in terms of the fields U and W: 

VUVW exp [-G [U, W] /k B T] 



with 



G [U, W] 
k B T 



<FapV ln Q A [{u- 
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•WO/2] <t>BpV ^ QbW- 
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(24) 



(25) 



(26) 



Alternatively, we could have started with the free energy functional (Q) and integrate out the Gaussian variables <&a 
and (f>s- So far no approximations have been used, but to find an energy functional that only depends upon W we 
now employ a saddle point approximation with respect to U: 



SG[U, W] 
5U 



= 



^(r)+^(r) = l 



(27) 



In 



Here we use the definition, see Eqs.(^) and (||): (j)* A = F° r ^b tne equivalent definition applies 

equilibrium the field variable U — Wa + Wb is conjugated to the overall density of the system, which is constant in 
an incompressible mixture. We therefore believe that the influence of this approximation on the description of the 
system through the field W is very small. We shall discuss this in detail below. 

If we replace U with the U* [W] that fulfils this constraint we end up with a free energy functional that only depends 
on the field variable W: 



G[W{v)\ =P Vx P 
k R T A N 



^ W 
v 4 x iV 



(FapV ln Q A [{U* 



-WO/2] 



N 



<t> B pV ^ Q B [(U*-W)/2] 



»A 



N 



n B 



(28) 



where we have used the fact that adding a constant field £ to U* does not change the value of G[W r (r)]. We chose £ 
in such a way that J v d 3 r U* = 0. 

A difficulty in describing the system with the order parameter W is the interpretation of the field fluctuations in 
terms of the physical density fluctuations. We can calculate the averages of the microscopic densities after the saddle 
point integration over U. These averages are marked by the subscript EPD. To this end, we introduce a local exchange 
potential A/i, which couples to the microscopic density difference 4>a — 4>b'- 
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with 



4> A pV ^ Q A [{U* + W + Afi)/2] _ BP V Q B [(U* -W-A(i)/2] 



4\N N n A N n B 

Thermodynamic averages of the microscopic density difference are obtained via functional derivatives: 

27V f SZ[A/j] 
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(30) 
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Similarly, we can calculate the fluctuations of the total density, which are induced by the saddle point approximation. 
1 
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Moments of the total density averaged over the field configurations of W are given by: 



(32) 
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(33) 



(34) 



These equations describe the actual fluctuations of the microscopic composition of the system after the saddle point 
approximation, i.e., of the EPD method. Having performed a saddle point integration in U , we have ignored fluc- 
tuations, and we cannot expect Eqs.||l]) or (34) to be accurate. Indeed, while the incompressibility constraint is 
fulfilled on average Eq. (Q) demonstrates that the saddle point approximation leads to spurious fluctuations of the 
total density. In Appendix p] w e use the random phase approximation (RPA) to evaluate Eq. ( |3l| ) and show explicitely 
the deviations between Eq.(pl|) and the well known RPA-structure factor. 

Alternatively, we can deduce the exact averages from the full free energy functional G[U, W] in Eq.(|2^) by intro- 
ducing a local exchange potential Afj, like in Eq.(E9|). After a variable substitution VF+A/u— > W this leads to: 



G[U,W,A(i]=G[U,W}--t I d 3 r 



-A11 2 + 2AfiW 
4%N 



(35) 



With this free energy functional we obtain the exact averages of the microscopic densities: 
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(j> A (r) -4> B (r))uw = — \- T (W) uw (36) 
X N 

—2i)(r — r') 1 

([^(r) - MrWA(r') ~ MOIW = { - 1 + -—^(W(r)W(r')) uw (37) 

px {x N r 

Recently, Ganesan and Fredrickson Q] have used a complex Langevin method to sample the fluctuations of both 
fields U and W and have obtained the average A-monomer density as (4>\)uw- Since U has to be complex to make 
the last term in Eq.(|24|) a proper representation of S($a + 4>b — 1) individual contributions to this average also have 
an imaginary part and the numerical procedure is quite involved. 

We expect the saddle point integration over U to be accurate whenever G[U, W] can be well approximated by a 
parabola in U —U* . In this case the fluctuations of W are only very little affected by the saddle point approximation in 
U and the fluctuations of W in the EPD method will closely mimic the fluctuations of W of the exact partition function 
(p4|). Hence, we can use Eq. (p7|) with (W(r)W(r')) uw (W(r)W(r')) (i.e., after the saddle point approximation) to 
obtain a very good approximation for the structure factor. In Appendix |^ we confirm that in RPA the fluctuations 
of W are, of course, not affected by the saddle point integration over U. Therefore we use the Fourier transform of 
Eqs.(^) and ( |3~7j ) in our calculations: 

(Mq)-Mq)) = -^(w(q)) (38) 

(\Ml) - Mq)\ 2 ) = + (-^2 (l^)l 2 ) (39) 

Because the composition is conserved and (4>a — 4>b) ~ we expect the order parameter W with which we 

are now describing our system to also be a conserved quantity. Therefore the dynamics of W are given through the 
relaxational dynamics of a model B system, referring to the classification introduced by Hohenberg and Halperin W4| . 



dW(r) 



= V r / A(r, r')V r ^w(r') + r?(r, t) (40) 
Jv 



dt 

with the chemical potential being the first derivative of the free energy with respect to the order parameter: 

^ (r) = 5 ^WW = J^ iw + xN [ ^ (r) - ^ (r)]) (41) 

The Fourier transform of this new diffusion equation is: 

= ~ Hq)q2 2]vbv [W{q) + xN{ ^ q) - ^ m + ^ q) (42) 

r) is white noise that obeys the fluctuation-dissipation theorem. The method using this diffusion equation we refer to 
as the EPD method 

We have found a diffusion equation that describes the dynamics in terms of the external field W = Wa — Wb and 
leads to the right physical equilibrium. A similar equation without noise has been derived by Maurits et al. pq| . The 
question to be asked is whether this dynamics represents any actual physical dynamics and how the choice of the 
Onsager coefficient influences the dynamics of the densities. It can be shown pq] , see also appendix that using local 
coupling in the EPD method is a good approximation for Rouse dynamics. The Onsager coefficient that we would 
have to use in DSCFT to reproduce Rouse dynamics is given in Eq.(|l^). The equivalent (local) kinetic coefficient in 
the EPD method for the same dynamics of the densities is: 

Aepd = -2 X ND (43) 



For the EPD calculations again the Fourier expansion of Eq.(|22|) is used. After having found the initial fields that 
create the given densities with the methods used for the DSCFT the chemical potential uw is calculated according 
to Eq.(ETt). [i w is then plugged into Eq.([42|) to find the time derivative of the difference in the fields. Thereafter 



15t~ = — ^15t~ ^ s integrated via the simplified Runge-Kutta method. After we have found the new W = Wa — Wb , 
we need to find the variable U* to make sure the incompressibility constraint 4>* A + (fi B = 1 given through the saddle 
point approximation ( p7| ) is fulfilled again using the Newton-Broyden method. The new fields lead us to the new 
chemical potential to calculate ^f- — d ^ t A — g t B and so forth. 

The method of the EPD has two main advantages compared to DSCFT: First of all it incorporates nonlocal coupling 
and secondly proves to be up to an order of magnitude computationally faster. There are two main reasons for this 
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huge speed up: In EPD the number of equations that have to be solved via the Newton-Broyden method to fulfil 
incompressibility is just the number of Fourier functions used. The number of equations in DSCFT that have to be 
solved to find the new fields after integrating the densities is twice as large. On the other hand, comparing diffusion 
equation (|l9| ) used in the DSCFT method with equation (|4p| ) in EPD it is easily seen, that the right hand side of the 
latter is a simple multiplication with the squared wave vector of the relevant mode, whereas the right hand side of 
Eq.(|l9|) is a complicated multiplication of three spatially dependent variables. 



III. MONTE CARLO SIMULATIONS 



A. Bond Fluctuation Model 

The Monte Carlo simulations presented in this study make use of the bond fluctuation model (BFM) (|32"| , j33| , which 
is a coarse grained lattice model, that incorporates the relevant features of polymers. These are connectivity of the 
monomers along a chain, excluded volume of the segments and thermal interaction between monomers. In this model 
each effective monomer occupies a cube of the lattice and blocks the eight sites at the cube corners for other monomers. 
Monomers of a chain are connected by one of 108 possible bond vectors of length 2, y/E, y/E, 3, or y/lO measured in 
units of the lattice spacing. (All following lengths are assumed to be given in these units unless an explicit unit is 
given.) These bond vectors are chosen to ensure the excluded volume condition, that makes sure that they do not cross 
each other during their movement. This large number of possible bond vectors allows 87 different bond angles, that 
provide a good approximation for continuous connectivity between the monomers of the chain. Each of the effective 
monomers represents 3-5 real chemical repeat units |45| , |46| . The number density of the occupied sites is chosen to be 
p = 1/16, which reproduces the properties of a polymer melt well. Interactions between the monomers are modelled 
through a square well potential with monomers of the same kind attracting and monomers of different kinds repelling 
each other. The interactions are chosen to be symmetric and to act inside a radius that extends over the first peak 
of the pair correlation function. This means the interactions act up to a distance of y/E which is equivalent to the 54 
neighbouring cubes of a monomer. 

ik B Te^0 tor y/E 
-e AA = -e BB = e AB = j Q {oi r > ^ ^ 

The moves used to simulate a purely diffusive movement of the monomers are local random monomer hopping 
moves, where one tries to move a randomly picked monomer to a neighbouring lattice site. 



B. Comparison between SCFT Calculations and Monte Carlo Simulations 

The Monte Carlo simulations were carried out with chains of 64 effective segments in a box with length L = 
160. With an overall number density of p = 1/16 there are 256 000 particles in the system. In the Monte Carlo 
simulations however the choice of N — 64, which is equivalent to a polymerisation of 200-300 in real polymers, is 
a compromise between the possibility to compare the simulation results with mean field results and the largest still 
sensibly manageable amount of computational resources, because an increase in the polymerisation leads to both an 
increase in the length scale and a slowing down in the kinetics of phase separation. To actually compare Monte Carlo 
simulations with SCFT calculations the parameters of both models have to mapped onto each other: 

In SCFT the only present length scale is the end-to-end distance R e of the polymers. It can be measured directly 
in the Monte Carlo simulations R e = y/Nb = 25.12. This gives us the length of the system to be L = 160 = 6.35i? e . 

In the SCFT calculations only the combination of Flory Huggins parameter and chain length enters. This sets 
the temperature scale. The Flory Huggins parameter \ can be calculated from the interaction parameter e of the 
square well potential defining the interaction between monomers in the Monte Carlo simulations with the following 
relation 

X = -r-~T^ z M e AB - \ {tAA + ess)] = 2z off e (45) 

KB 1 2 

z e ff is the effective coordination number in the bulk, i.e. the average number of intermolecular contacts per monomer. 
We hereby speak of contacts, if the distance between the monomers is smaller than y/E. 

The average composition <f> A is a parameter of the SCFT, the total number density p of monomers is only required 
if fluctuations are considered. 







The single chain diffusion constant D, which can be extracted from the Monte Carlo simulations by measuring the 
mean square displacements of the chains, gives the time scale r = R\jD — 1.5 x 10 7 Monte Carlo steps (MCS). r is 
constant, because D and R e are almost independent of time and composition Q. 

The length of the chains used in the Monte Carlo simulations is N — 64. This is somewhat larger than the 
entanglement length N e 32 ^5j,^(|, therefore, we are in a cross over regime between Rouse dynamics and reptation 
1(491, j50[ . This means, when we are comparing dynamic mean field results with Monte Carlo results we expect to find 
a reasonable agreement if we regard Rouse dynamics. 

This identifies all parameters of the SCFT calculations (without noise). If we neglect fluctuations, systems with 
identical x^V, R e , and composition but different degree of interdigitation pR^/N |34|,[5l| and statistical segment length 
b = R-l/N give identical results. The degree of interdigitation controls the strength of fluctuations and mean field 
theory is believed to be correct in the limit pR^/N — > 00. The statistical segment length sets the smallest length 
scale the Gaussian description of polymers is valid. If we were interested in the structure on smaller length scales we 
would have to use a different chain model (e.g., worm-like chain |52,p3]). 




In the following sections we regard the early stages of spinodal decomposition after a quench from the one phase 
region with xN = 0.314 into the miscibility gap with xN — 5 for a symmetric binary polymer mixture. For much larger 
incompatibility the width of the interface becomes of the order of the statistical segment length, and properties on this 
length scale cannot be described by the Gaussian chain model. For smaller incompatibilities -in the vicinity of the 
critical point xN = 2- composition fluctuations are very strong (i.e., non-Gaussian) and the mean field approximation 
becomes worse. We will first show some general results for the dynamical mean field theory. We will then carry on 
to compare these results with Monte Carlo simulations showing what role the Onsager coefficient plays. After having 
neglected random statistical fluctuations in SCFT so far we will explain how fluctuations are implemented into the 
two dynamic mean field theories and discuss their influence on the dynamics. 



If a system is quenched from the one into the two phase region the linearized theory of spinodal decomposition 
P-|TT| predicts that fluctuations with wave lengths A larger than a critical value A c , i.e wave vectors q below a critical 
value q c , start growing spontaneously. This is illustrated in Fig.[l] showing results obtained through DSCFT with local 
coupling in a one dimensional system. Before the quench at time t = we have a homogeneous mixture with random 
statistical density fluctuations. After the quench these fluctuations are amplified. It is apparent that the resulting 
density profile is well describable with plane waves, so that the set of basis functions we have chosen proves to be not 
only technically convenient for the calculations but also describes the physical phenomena well. At later times we 
see that the mode with the wavelength about a third of the system size is amplified most until at even later stages a 
saturation inside domains takes place until these domains are separated by sharp interfaces. Later stages of demixing 
are not appropriately described within our purely diffusive model because hydrodynamic mechanisms and random 
fluctuations are neglected but play a dominant role. 

Linearised Cahn-Hilliard-Cook theory P~pl| predicts an exponential behaviour for the density modes 4>a (q) ~ 
e R(q)t as i on g as th e difference 4>a (r) — 4>a between the actual density and the average density is small. Starting 
with a homogeneous mixture before the quench this requirement is fulfilled during early demixing. Our density 
coefficients </u(q) are equivalent to the actual, but discretised, density modes. In Fig.g some of the coefficients cf>A (q) 
are displayed versus time. These results were obtained for a three dimensional system using DSCFT with a local 
Onsager coefficient. The exponential behaviour is obviously well reproduced for early times. For larger values of 
q the exponential behaviour changes earlier. This is also in agreement with experimental results p4j| . The faster 
growth of density modes with smaller wave vectors leads to the creation of small domains on a length scale ~ 1/q. 
These domains then cause density modes with wave lengths smaller than the extension of the domain to be damped. 
Fig.^ shows the corresponding relaxation rate R(q) versus q. We see that for wave vectors with a value below q c an 
exponential growth of the density modes sets in whereas for q > q c initial fluctuations are damped. In the growth 
region there is a maximum. The two dashed lines framing our results are the relaxation rates for temperatures just 
below the critical temperature (weak segregation limit WSL) and for very low temperatures (strong segregation limit 
SSL) as they are given by the Cahn-Hilliard-Cook theory: 




IV. RESULTS 



A. General aspects of spinodal decomposition 




(46) 
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with 



v2fc — 

li (47) 

= ^ {x N/2-l)^ 

The last line is only valid for symmetric mixtures with (f>A = 4>b — 0.5 when the spinodal point lies at xsN = 2. k is 
found to be 18 in the WSL j23| and 12 in the SSL |p5|| . The maximum growth rate is found for q m — q c /\/2. With 
\N = 5 lying between the SSL and the WSL our result was to be expected. The influence of Rouse dynamics on the 
relaxation rate is displayed in Fig.||. Here we have presented results for a two dimensional system using DSCFT with 
Onsager coefficients describing both local kinetics and Rouse kinetics for an approximately homogeneous mixture and 
EPD. The EPD and the DSCFT results using Rouse dynamics are virtually identical, but remembering that the EPD 
method is computationally much more favourable, this method should be used when Rouse kinetics is considered. 
Local and nonlocal coupling lead to distinct differences in the relaxation rates: in the range < q < q c nonlocal 
coupling leads to a reduction of the relaxation rate and a shift of the maximum to smaller values of the wave vector. 
Also in the region q > q c the initial density fluctuations are not as strongly damped as for local dynamics. 

Note that the difference in growth rate between the local dynamics in the WSL and SSL is very similar in magnitude 
to the difference in growth rate between the full SCFT calculation using a local or a Rouse Onsager coefficient. 
This demonstrates that the square gradient expression for the free energy is not sufficiently accurate: quantitative 
deviations from experiments or simulations might be either due to the additional approximations of the square gradient 
approach or the wave vector dependence of the Onsager coefficient. To compare simulations and theory quantitatively 
and to extract the interplay between single chain dynamics and the kinetics of collective composition fluctuations a 
quantitative self consistent field calculation is required. 



B. Comparison between Monte Carlo results and dynamic mean field results 

The Monte Carlo simulations were performed on a Cray T3E using a trivial parallelisation scheme running 64 
configurations in parallel to achieve good statistics. 5 400 000 MCS were performed which is equivalent to 45 days of 
CPU time per processor. The EPD and DSCFT calculations, the results of which are presented in this section, are for 
a three dimensional system using 7x7x7 = 343 basis functions. The equivalent to 12 700 000 MCS were performed 
taking approximately 65 days on a Cray J90. 

The global structure factor defined by 



S(q,t) = ( / d 6 r [«^(r,i)-MM)]e iqr > ( 48 ) 



v 



is an important experimentally measurable quantity for the description of the phase separation process. In Fig.^ the 
global structure factor is plotted versus wave vector q for different times, thin lines denoting SCFT results, wide lines 
with the same symbols corresponding Monte Carlo simulation results. (Because the SCFT calculations could only be 
performed for a single starting configuration instead of 64 like in the Monte Carlo simulations, for the initial time 
the global structure factor of a homogeneous mixture given by RPA |^] was used. Global structure factors for later 
times were extracted through the exponential time dependence of the density modes in the SCFT calculations.) Part 
a) compares Monte Carlo with local kinetics and part b) Monte Carlo with Rouse kinetics. Fig.||a) clearly shows no 
quantitative agreement between DSCFT with local kinetics and the Monte Carlo simulations: The peak in SCFT 
grows much quicker and also the position of the peak is too far right. As has been mentioned before, including Rouse 
dynamics leads to a reduction of the relaxation rate in the growth region and the wave vector with the maximum 
growth rate is shifted towards a smaller value. Fig||b) satisfies these expectations: The position of the peaks for SCFT 
and Monte Carlo almost coincide and the growth rates of the peak are much closer to each other although SCFT 
still overestimates the growth rate. A more detailed comparison is possible if we plot the corresponding relaxation 
rates versus wave vector. As we have seen before the modes of the density calculated with the dynamic mean field 
theory follow a clear exponential behaviour. To derive the relaxation rates of the global structure factor, we need 
to plot the modes of the global structure factor on a logarithmic scale versus time. This is done in Fig.j| for some 
randomly chosen values of q. The modes with a q value smaller than the critical value q c derived through DSCFT 
show a clear exponential behaviour. This behaviour changes with time, especially the bigger the value of q of the 
growing mode the sooner the change in the exponential behaviour takes place. This is qualitatively in agreement 



with the results we obtained through dynamic mean field theory, see section [V A. Modes with q > q c are more or 
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less constant in our figure. Therefore an estimation of q c from the Monte Carlo results alone would suffer from some 
ambiguity: of course, this is not only because the accuracy of the Monte Carlo results is limited, but is a matter 
of principle § . The interplay of fluctuations and nonlinear effects has the consequence that a well-defined q c does 
not exist. For times earlier than the displayed interval random fluctuations influence the behaviour of the modes so 
strongly, that an exponential behaviour is not visible. The relaxation rates resulting from the fits to the points in 
Fig.^j for the two different time intervals indicated through the solid lines and the corresponding mean field relaxation 
rates are displayed in Fig.]?]. Only the exponential behaviour of the earlier of the two marked time intervals may 
be interpreted as the expected behaviour of early spinodal decomposition. For later times a gradual change away 
from the exponential behaviour sets in. The second fit shows an apparent exponential behaviour, because the time 
interval is too small to resolve the change. Hence the resulting 'relaxation rate' for the later time interval may only 
be treated as an indication for the deviation from the earlier exponential behaviour. Part a) shows the relaxation rate 
for DSCFT with local dynamics in the time interval < t < 1 x 10 7 MCS and the Monte Carlo results for the two 
time intervals 7.25 x 10 5 < t < 1.5 x 10 6 MCS and 1.92 x 10 6 < t < 2.63 x 10 6 MCS; part b) shows the same results 
but regarding Rouse dynamics through the use of EPD and DSCFT with Rouse dynamics. The inset figures are the 
corresponding Cahn plots displaying R(q)/q 2 versus q 2 . As is obvious from the behaviour of the global structure 
factor, local dynamics gives a relaxation rate that is much too large and the wave vector of the maximal growth rate 
is also too big. The Cahn plot of the mean field relaxation rate shows the linear behaviour, which is expected for 
local coupling. This is in strong disagreement with the Monte Carlo results. For Rouse dynamics we actually find an 
almost quantitative agreement for earlier times in the region of the positive growth rate. For later times, however, 
the relaxation rate decreases in the Monte Carlo simulations - this has also been observed in experimental studies 
p7| -, while SCFT still shows the same exponential behaviour as is found for earlier times. In the Cahn plot we see 
that Rouse dynamics leads to a nonlinear run of R(q)/q 2 versus q 2 . This is obviously also the case for the simulation 
results and is in agreement with earlier simulations Jl5|-p^7t and experimental observations |27]-|2£]. This nonlinear 
behaviour of the Cahn plot is related to the fact that we consider a deep quench, for which q c R e > 1 [Z2|. For a shallow 
quench, for which \N exceeds the critical value XcN — 2 Q only slightly, one has q c R e < 1 [§2j, and then the theory 
would yield a linear Cahn plot, also consistent with corresponding observations ||. The latter case is less interesting, 
however, because then the polymer mixture is to a large extent equivalent in behaviour to a fluid mixture of small 
molecules, and there is no longer an effect of internal Rouse relaxation modes on the phase separation dynamics in 
this limit. 

If we now compare our results for larger wave vectors we see big discrepancies independent of the chosen dynamics 
in SCFT: At q c the structure factor is independent of time in the early stages of the SCFT calculations, i.e. structure 
factors at different times exhibit a common crossing point at q c , see the inset graphs of Fig.||. No such intersection 
occurs in the simulation data. Mean field theory damps the density modes with q > q c while Monte Carlo simulations 
lead to a relaxation rate fluctuating around zero. This behaviour is also seen in the global structure factor, see FigJ^: 
The right side of the SCFT peak decays fast, while the Monte Carlo peak is much broader with a slower decay. Both 
the earlier change of the exponential behaviour and the form of the relaxation rate for larger wave vectors are the 
result of random fluctuations. As was expected, the influence of fluctuations during the very early stages of spinodal 
decomposition on the growth of the density modes is rather small, but they are crucial for smaller wave lengths 
and determine the change in the exponential behaviour, because random fluctuations cause some modes to reach an 
amplitude where the non-linear regime sets in early. 



C. The influence of random fluctuations 



The diffusion equations we have used so far are completely deterministic, but obviously in all dynamic processes 
random statistical fluctuations are present. To regard these fluctuations we have to add a random force i] to our 
diffusion equations (19) and (Eoh that is linked to the Onsager coefficient through the fluctuation-dissipation theorem: 



(V(r)) = 

(rj(r, t)r]{r', t')) = 2fc B TA(r, r')W 2 S(r - r') S(t - t 1 ) 



(49) 



or Fourier transformed: 



(%) = o 

Mt)V- q (t')) = (kl 2 ) - 2k B TA(q)q 2 S(t t 1 ) 



(50) 
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In our calculations the diffusion equation is integrated through discrete time intervals St that are determined through 
the Runge-Kutta scheme. In DSCFT, for example, we use the Langevin equation: 



4>A, q (t + St) = (f> A ,q(t) + 6tA{q)q 2 n Q + f q (St) (51) 
with fq(St) expressing random fluctuations that obey the fluctuation dissipation theorem [p^prl : 

f q (6t) = ^/2A(q)q^V6ir (52) 

r is a random number with the properties (r) = 0, (r 2 ) = 1. 

The analogue f q (St) is used to perform EPD calculations with fluctuations. The difficulty in this case is to interpret 
the resulti ng fl uctuations of the external field W = Wa~Wb in terms of the density $ = {4>a~ §b)J2l- As we have seen 
in section II C the field and the density fluctuations are linked to each other through equations (Bq) and (p9[ ). 

To ensure that the way the fluctuations are included is correct and the validity of these Eqs.fl38|) and fl39|) is given, 
we consider a homogeneous, %iV = 1.8, one dimensional system with length L x = 6.35i? e using both DSCFT and 
EPD. From the DSCFT calculations we derive (l^l 2 ) for each wave vector while the EPD method leads to (lIFql 2 ). 
In both cases we averaged over 10000 snapshots of the density or the field made during their time evolution. In Fig.|8| 
(|$ g | 2 ) and (iVFgl 2 ) are plotted versus q. The solid line is the expected result from RPA. For larger values of q we 
find good agreement, but the amount of used configurations was not enough to find reliable results for smaller values 
of q because the correlation times are much longer due to the diffusive dynamics. If we plot (|$ g | 2 ) versus ( | Vt 7 ^ | 2 ) , as 
is done in Fig.^, we find good agreement of our data with the expected linear behaviour. 

After having proven that the above treatment of random fluctuations leads to the expected behaviour we focus on the 
influence of fluctuations on spinodal decomposition. Monte Carlo results show that the exponential behaviour changes 
earlier than in SCFT calculations without fluctuations resulting in a reduced relaxation rate. The relaxation rates for 
different time intervals obtained through the Monte Carlo simulations and the EPD method are presented in Figjicj. 
The EPD results were obtained by averaging the time evolution of the fields of 64 two dimensional configurations. For 
earlier time intervals the relaxation rate is quantitatively very similar to the rate without fluctuations in the region 
below q c . For later times however when calculations without fluctuations still show the same exponential behaviour 
the relaxation rate is reduced as is also seen in the Monte Carlo simulations. In the range above q c the modes of the 
fields or densities do not follow an exponential behaviour but fluctuate around zero. Consequently the relaxation rate 
is not well defined in this region leading to strong fluctuations of the relaxation rate in Figjl^. For very early times 
when the density modes are still of the order of the fluctuations of the homogeneous system it is also not possible to 
see an exponential behaviour because the density changes caused through random fluctuations conceal the growth of 
the modes during very early spinodal decomposition. 



V. SUMMARY 



In this study we analysed the influence of single chain dynamics on the collective diffusion during early stages of 
spinodal decomposition in a symmetric binary polymer blend. We used the SCFT for polymer mixtures to explore 
two versions of dynamical mean field theory. The single chain dynamics enters these descriptions through an Onsager 
coefficient. The first method we call DSCFT propagates the densities in time and gives us the possibility to model 
both local dynamics and approximately the nonlocal dynamics we expect for the Rouse model. In DSCFT the 
correct treatment of Rouse dynamics would involve the calculation of the pair-correlation function at every time step, 
which is a computationally rather expensive task. On the other hand during early stages of demixing the mixture 
is only weakly inhomogeneous so that the use of the pair-correlation function of a homogeneous mixture, which is 
analytically known, serves as a sufficient approximation. In the second method the instantaneous configuration is not 
described by the densities but by the effective external field (EPD). We find a Langevin equation for the external 
field which using 'local' kinetics in the fields is found to describe a polymer mixture with Rouse dynamics. Apart 
from 'automatically' including Rouse dynamics this EPD method has the big advantage of being up to an order of 
magnitude computationally faster than the DSCFT method. 

First numerical calculations with these methods for a quench from the one phase, %iV = 0.314, to the two phase, 
xN — 5, region neglecting random fluctuations show a clear exponential behaviour of the density modes, as was 
expected for a mean field description. The relaxation rate of the density mode is strongly influenced by the choice of 
the Onsager coefficient: In the growth region Rouse dynamics reduces the relaxation rate compared to local dynamics 
and the position of the maximum growth rate is shifted to smaller values of q. For q > q c however Rouse dynamics 
causes the modes to be damped less quickly. We also find good agreement for early stages between the DSCFT 
method using the pair-correlation function for the homogeneous system and the EPD method. 
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To quantitatively test these mean field predictions we compare them with results obtained through Monte Carlo 
simulations employing the bond fluctuation model. The chains used in these simulations are expected to show Rouse 
behaviour because of the chosen chain length N = 64. The comparison is possible without any adjustable parameter. 
We compare the global structure factor, which is the experimentally accessible quantity, and the relaxation rates. 
Local dynamics in DSCFT overestimates the growth of the global structure factor by far, but the agreement is better 
for Rouse dynamics especially for earlier times. Neglecting random fluctuations in our mean field calculations proves 
to be justified for wave vectors with positive relaxation rates and earlier times, but should be included to investigate 
later times. 

Fluctuations can easily be included in DSCFT and EPD. The difficulty in the EPD method is the fact that the 
field fluctuations have to be interpreted in terms of the physical density fluctuations. We can find a relation between 
the field and the density fluctuations. EPD calculations with fluctuations lead to an earlier change in the exponential 
behaviour of the density modes as was also the case for the Monte Carlo simulations. The missing dampening of the 
modes with q > q c as is found in the simulations is also reproduced. 

We have seen that the single chain dynamics has a pronounced influence on the collective dynamics of a polymer 
mixture. Comparing quantitatively Monte Carlo simulations and dynamical mean field theory we have validated 
the mean field calculations. Note, however, that we have considered a deep quench far below the critical point; for 
shallow quenches close to the critical point mean field theories are not expected to be accurate. The later stages 
of spinodal decomposition are not accessible with either method. During later times hydrodynamical interactions 
become important. Lattice model Monte Carlo simulations lack a hydrodynamical mechanism. In dynamical mean 
field calculations hydrodynamic coupling can be included |^8|,[5^], and, hence, they can be extended beyond the validity 
of the lattice model. 
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APPENDIX A: REPRESENTATION OF ROUSE DYNAMICS THROUGH EPD 



As mentioned before it is possible to show that the EPD method using local coupling is a good approximation for 
reproducing Rouse dynamics of the physical densities. The derivation of the EPD method we present in this section 
was introduced by N.M. Maurits et al. and can be found in reference p6| . 

Following the method we used to include dynamics in the SCFT to achieve the DSCFT method we employ a saddle 
point approximation in the external fields which leads to a bijective relation between the external fields Wa, u>b and 
the densities 4>a, i>B- This means we can choose with which of the two variable sets we would like to calculate. On 
the other hand the pair-correlation function which is part of the Rouse Onsager coef ficient can be expressed as the 



functional derivative of the density with respect to the external potential, see also Eq.(B13) 

(A1) 

To calculate in w space we have to transform the time derivative of the densities according to the chain rule: 

d A^A = f HA(r,t) d WA (r',t) = Ja f , ^ A (r',t) 3 

dt J Sw A (r',t) dt N J oy ' ' dt y ' 

Combining this equation with the diffusion equation for Rouse dynamics (|l8| ) this leads us to 

"I / ^o(r,r')^f^^ = ^V r ^P (r,r0V r ^(O^ (A3) 

Using the approximation 

V r P (r,r') = -V r ,P (r,r') (A4) 
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one easily arrives at 



- / Po(r,r') 
Jv 



„d WA (r> lt ) drl3 = D 



leading to an equation of motion for the external fields 

dw A (r) 
dt 



[ P (r,r')V^ A (r')dr' 3 
Jv 



-DV 2 p A (v) 



(A5) 



(A6) 



Approximation (A4) is obviously exactly valid for a homogeneous mixture, because the pair-correlation function only 
depends on the distance |r — r'| between two points. In the inhomogeneous case, if r and r' are in different phases and 
neither of them in the interface, V r Po(r, r') and V r 'Po(r, r') are of different sign making this approximation justifiable 
even if the actual values differ. 

Because of the incompressibility constraint p A and ps are not independent of each other. There is only one 
independent chemical potential p = p A — pg so the equation of motion for the external fields which is to be used in 
our case has the form: 



dw A (r) 9ws(r) 



dt 



i)t 



-DV 2 p(r) 



(A7) 



If we compare this equation with equation ( }40| ) we see that densities evolving in time according to Rouse dynamics 
are well described through the EPD method if a local kinetic coefficient as given in Eq.(0) is used. 



APPENDIX B: RANDOM PHASE APPROXIMATION FOR THE FLUCTUATIONS IN EPD 



The single chain partition function Q is defined through (compare with Eq 



(Bl) 



where </>i(r) = N/p J ds S(r — r(s)) denotes the single chain density. We now expect the system to be only weakly 
inhomogeneous, i.e., meaning the density and the external field only differ a little from the average value: 0i(r) = 
N/pV + 5(f>i(r); W(r) = W + SW(r). The density and the external field are presented as a Fourier expansion: 



k(r) = N/ P V + J2ke' 

q#0 



qr 



V 



d 3 r 0i(r)e 



(B2) 



This Fourier expansion is now inserted in Eq.(Bl) 



Q= / X>[n]7>i[r]exp 



■■exp[-WA] J ©[riFiM jl-^ 
: exp [-W A ] QoU-jfY, W ^-i 



q#0 



q,q'#o 



(B3) 



p 2 V 2 



2N 2 ^ 



W q W q 4-q< 



Qo denotes the partition function of a single chain without an external field. The average in the last line is to be taken 
over all chain configurations that are possible when there is no external field present. Because the average deviation of 
the density (5<ft)o from the average value is zero if there is no external field, (0 q )o = is valid. The average (</>_ q 0_ q <)o 



is given through the single chain structure factor S'o(q) of a Gaussian chain: ((/>. 
Neglecting higher terms, we obtain the RPA-result for the single chain partition function: 



— N 
-q'/o — TSyir 



S (q)<5- 



Q RPA = exp [-W A ] Ooexp 



^5>(q)|wg 



q#0 



q,-q'- 



(B4) 
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Obviously a corresponding expression is valid for a B polymer. These RPA single chain partition functions are plugged 
into Eq.fPq) leading to: 



G[U,W] _ tApV^QAfl tBpV ^QBfi | pVx [ (<Fa-4>b) P VW | pV A - r2 



N 



tia N n B 
P V 



q#0 q#0 



Ei^ 



2A 

(Wq + ^q) 



4iVxiV 

_ <?WU 

2N 2 



E^(q) 



q#0 



(B5) 



Regarding only the wave vector dependent parts of the free energy we find: 



G [U, W] 



pV 
2N 2 



pV 
2N 2 



E 

q^O 

pV_ 

2N 2 

E 

q#0 



So(q) 1 
4 2 X 



5o(q) 



If/a 



E 



q#0 



-5 (q) (C7 q W- q ) + G hom = 



go(q) 1 1 
4 2 X 4 l 



B ) 2 S (q) l^ql 2 + 



(B6) 



pV f^o(q) 



2iV 2 



U7 q + (</»A-0s)W q f + G 



q#0 



T hom 



We use this free energy to evaluate the partition function of Eq.fl2q). Following the procedure we used before, we 
employ a saddle point approximation with respect to the field U: 



SF[U,W] 

su 







u* 



(B7) 



Since the free energy is quadratic in the RPA approximation the saddle point integration is equivalent to the functional 
integration over U. This again leads us to an expression for the free energy G[W] only depending on the external 
field variable W: 



G[W] _ pV 



k B T 2N 2 
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q#0 



jSo(q) 
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2X 
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(B8) 



Now the average ( | VFq | 2 } can be calculated: 



<|wg 2 ) = 



N 



2X 



pV l-2x0 A &?So(q) 



Using our result (p 



(I^Aq 



1 



, 2 <W> 



(B9) 



(BIO) 



pV X (XN) 2 

we express the density fluctuations in terms of the field fluctuations and recover the well-known RPA expression W2[ 



Aq - <PBq| 



P v 
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Using the RPA single chain partition function (B4), we calculate 

t> A V SQ A 



= - 

ft to 



Qa SW A (r) 



= *A-^£So(q)Wk q e* 



q/0 
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6W A (r') " NV 
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jq(r-r') 



(Bll) 

(B12) 
(B13) 



q#0 



1G 



The last equation is equivalent to Eq. Al but due to the particle conservation the q = contribution has to be taken 
out of the sum. <\>a (q = 0) is just the average overall density <j>A that cannot change if the external field is altered. 
With these expressions we obtain for the "literal" fluctuations in the EPD method according to Eq. ([Il]) : 



N f i 
(|<M q " 0Bq| a )BPD = (\4>m ~ 0B q | 2 > - ^72 y d 3 r' e^*"' 

8 x S 2 (q)<^| 5 (q) 



5W A {r') 5W B (r') 



pV(l-2x<f>A<f>BS (q)) pV 
= « + ^f^ CB14) 

Generally, the deviation for the RPA result is of similar magnitude as the RPA structure factor itself. For a symmetric 
quench (j>A = \, however, we accidentally recover the RPA result. This example also illustrates that one can obtain 
the average of the composition by sampling the average of <j) A — <p%, i.e., the densities of single chains in the field 
configuration W, but one should not use this to calculate fluctuations. 
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FIG. 1. Density profiles at different times during demixing after a quench from \N = 0.314 to X-^ — 5 in a one dimensional 
system using 12 eigenfunctions. The upper panel shows the early stages when concentration fluctuations with wavelengths 
between a third and a quarter of the system size are amplified. In the lower panel later stages are displayed - the concentration 
inside a domain slowly saturates leading to sharp interfaces between the coexisting phases. Note the change of scale on the 
composition axis between the two panels. 
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FIG. 2. Several density modes displayed on a logarithmic scale versus time. The displayed values of q are given in units of 
1/Re- The results were obtained through DSCFT calculations in a three dimensional system of length L x =L y = L z =6.35i? e 
using 7x7x7 functions for a quench from \N = 0.314 to \N = 5. The expected exponential behaviour during early stages of 
demixing is well reproduced. 
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FIG. 3. Corresponding relaxation rate to the density modes displayed in Fig|| Below a critical wave vector q c the density 
modes are increased spontaneously. Modes with larger wave vectors are damped. As expected the results are found to be 
between the two limits of weak (WSL) and strong segregation (SSL) as given by the Cahn-Hilliard-Cook theory. 
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FIG. 4. Relaxation rates obtained through two dimensional DSCFT calculations using local and Rouse dynamics and EPD 
calculations. The DSCFT results using the pair-correlation function of a homogeneous melt and the EPD results are in good 
agreement. 




FIG. 5. Global structure factor versus wave vector for different times. Broader lines represent Monte Carlo results, thin lines 
with the same symbols the corresponding DSCFT results. Panel a) compares DSCFT with local coupling with the Monte Carlo 
simulations. Local dynamics obviously overestimates the growth rate and shifts the maximum growth rate to larger values. 
Panel b) compares Rouses dynamics with Monte Carlo results showing better agreement. The inset graphs show the behaviour 
of the global structure factor in the area of q c : While the mean field results lead to a common intersection point defining q c , 
the Monte Carlo lines do not cross each other in a single point making the definition of q c impossible. 
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FIG. 6. Global structure factor versus time for a few randomly chosen values of q. The values of the wave vectors are given 
in units of 1/R e . q c corresponds to the critical wave vector extracted from DSCFT. An exponential growth of the modes is 
found. Modes with a smaller growth rate change their exponential behaviour earlier than those growing faster. 




FIG. 7. Corresponding relaxation rates to Fig.g. Panel a) compares the Monte Carlo relaxation rates with DSCFT calcula- 
tions with local dynamics. Panel b) compares Monte Carlo results with EPD and DSCFT calculations with Rouse dynamics. 
For earlier times good agreement in the growth region is found but Monte Carlo simulations show an earlier change in the 
exponential behaviour. The inset figures are the corresponding Calm plots, where one displays R(q)/q 2 versus q 2 . 
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FIG. 8. The average (|$ ? | 2 > derived with DSCFT and 



| VKj | 2 ) derived with EPD displayed versus q. Both calculations are 
valid for a homogeneous, \N = 1.8, one dimensional system with length L — 6.35i? e . For larger values of q good agreement 
with the RPA averages is found. For smaller values of q too few independent configurations of the system were taken into 
account. 
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FIG. 9. The linear behaviour expressed in Eq.(^) is well reproduced, as we can see, when (I*,! 2 ), obtained through DSCFT, 
is displayed versus (IWql 2 ), obtained through EPD. The points to the right correspond to small wave vectors, those on the left 
to large wave vectors. The solid line is the corresponding RPA result. 
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FIG. 10. Relaxation rates obtained through Monte Carlo simulations and EPD calculations in two dimensions with random 
fluctuations for different time intervals. Fluctuations lead in both methods to an earlier change in the exponential behaviour 
of the increasing modes. 
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